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ABSTRACT 


We analyze spectroscopic and photometric data to determine the projected inclinations of stars in 
11 open clusters, placing constraints on the spin-axis distributions of six clusters. We combine these 
results with four additional clusters studied by Healy & McCullough (2020) and Healy et al. (2021) to 
perform an ensemble analysis of their spins. We find that eight out of ten constrained clusters (8096) 
have spin-axis orientations consistent with isotropy, and we establish a lower limit of four out of ten 
(4096) isotropic clusters at 7596 confidence, assuming no correlation of spins between clusters. We also 
identify two clusters whose spin-axis distributions can be better described by a model consisting of an 
aligned fraction of stars combined with an isotropic distribution. However, the inclination values of 
these stars may be influenced by systematic error, and the small number of stars modeled as aligned in 
these two clusters precludes the interpretation that their stellar subsets are physically aligned. Overall, 
no cluster displays an unambiguous signature of spin alignment, and 9796 of the stars in our sample 
are consistent with isotropic orientations in their respective clusters. Our results offer support for the 
dominance of turbulence over ordered rotation in clumps and do not suggest alignment of rotation 


axes and magnetic fields in protostars. 
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1. INTRODUCTION 


'The process of star formation begins primarily within 
giant molecular clouds (GMCs). Within GMCs, 
“clumps” of higher density substructure form with com- 
parable mass and length scales to those of open clusters 
(Shu et al. 1987, and references therein). Large-scale 
supersonic turbulence is thought to help shape these 
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clumps, which become gravitationally unstable and be- 
gin to collapse as they decouple from the overall turbu- 
lent flow in their GMC (Lada & Lada 2003, and refer- 
ences therein). Within clumps, material fragments into 
protostellar cores, which can accrete additional molec- 
ular gas before clearing the star-forming region and 
emerging as an open cluster. The spins of the result- 
ing stars are drawn from the angular momentum of the 
material from which they formed (McKee & Ostriker 
2007, and references therein). 

Within an open cluster, determining the orientation of 
stellar spin axes facilitates an inference of the balance of 
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turbulent and rotational energy that existed during the 
formation of the cluster’s stars. While the dominance 
of turbulent kinetic energy in a clump randomizes the 
resulting stellar spins, numerical simulations suggest a 
comparable contribution from ordered rotation could re- 
sult in spin alignment for stars > 0.7 Mo via protostellar 
accretion of material sharing the clump's global angular 
momentum (Corsaro et al. 2017; Rey-Raposo & Read 
2018). Previous determinations of stellar spin-axis ori- 
entations in open clusters have yielded mixed results, 
finding isotropic spins in the Alpha Per, Pleiades, and 
Praesepe clusters (Jackson & Jeffries 2010; Jackson et al. 
2018, 2019) and broadly aligned spins in NGC 6791, 
NGC 6811, (Corsaro et al. 2017) and Praesepe (Kovacs 
2018). 

Our past contributions are also consistent with mul- 
tiple explanations. Healy & McCullough (2020, Paper 
I) determined the projected inclination (sin?) values of 
NGC 2516 cluster members, finding that their distri- 
bution is consistent with either isotropy or moderate 
alignment. Healy et al. (2021, Paper II) performed a 
similar analysis of the Pleiades, Praesepe and M35 clus- 
ters. The results included the possibility of a subset of 
aligned stars in M35, but not at a statistically significant 
level. As seen with NGC 2516 in Paper I, Paper II de- 
termined a similarly degenerate scenario for the Pleiades 
and Praesepe, with the data supporting either isotropy 
or moderate alignment. 

In this paper, we broaden our study of stellar spins 
with an analysis of 11 open clusters, 10 of whose stars' 
Spin-axis distributions have not been determined pre- 
viously. Our analysis is again driven by the spectro- 
photometric method of determining projected inclina- 
tions, described by the following equation: 


vsini = E sin i (1) 
Above, v sin? is a star's projected rotation velocity, R is 
the radius, and P is the equatorial rotation period. 

In expanding our sample, we not only aim to test the 
longstanding assumption of isotropic spins and learn the 
characteristics of each individual cluster, but we also 
seek to synthesize the results from the 11 clusters in this 
work and four from our previous papers to challenge the 
degenerate solutions determined in our previous works 
using insights from the overall ensemble of clusters. 

With these goals in mind, the structure of this paper 
is as follows: Section 2 lists the sources of data we use 
to perform our analysis, which is detailed in Section 3. 
Section 4 presents our results for this expanded study. 
In Section 5, we interpret each cluster's results and dis- 
cuss the broader inferences about star formation that 


we are afforded by our series of papers. Finally, we offer 
concluding thoughts and future outlooks in Section 6. 


2. DATA SOURCES 


As in Papers I and II, we used the cluster membership 
tables of Cantat-Gaudin et al. (2018, CG18) updated 
with Gaia EDR3 data (Gaia Collaboration et al. 2020) 
to establish our initial lists of stars to study. Recent re- 
search has revealed large tidal tails extending from the 
central regions of open clusters (Kounkel & Covey 2019; 
Meingast et al. 2021). We continue to rely on the CG18 
catalog since the cluster membership probability for in- 
dividual stars in tidal tail regions is more uncertain than 
for stars in cluster centers, and past spectroscopic sur- 
veys did not targets regions outside the central regions 
to provide the necessary data for our study. 

We constructed spectral energy distributions (SEDs) 
beginning with queries of the following catalogs also 
used in Paper II: (CatWISE2020 (Marocco et al. 2020), 
2MASS (Skrutskie et al. 2006), Hipparcos/Tycho-2 
(Perryman et al. 1997; Hog et al. 2000) and GALEX 
(Martin et al. 2005)). Using Gaia match tables to find 
corresponding entries, we also added SkyMapper DR2 
(Onken et al. 2019) and SDSS DR13 (Albareti et al. 
2017) colors along with Pan-STARRS g-band photome- 
try (Chambers et al. 2016; Flewelling et al. 2020) to the 
SEDs when available. We also added Gaia G magni- 
tudes to each SED, with mean uncertainties calculated 
using the flux values and errors and added in quadrature 
to the uncertainty in the G magnitude zero-point (Evans 
et al. 2018). The quality flags we applied to provide re- 
liable magnitudes for fitting are listed in the Appendix 
of Hamer & Schlaufman (2022). 

All of the time-series photometry in our analysis is 
based on TESS images and products (Ricker et al. 
2015). Table 1 lists each cluster in this study along with 
its galactic latitude and longitude, distance via inver- 
sion of Gaia EDR3 parallaxes, age determined by this 
work, mass estimate when available, and the source of 
vsin? measurements. This table also provides the v sin? 
threshold we established for each dataset (see Section 
3.4). For SDSS DR17 sources, we utilize data from the 
APOGEE-2 spectrograph (Majewski et al. 2017). Data 
for Pozzo 1 and NGC 2547 originally come from the 
Gaia-ESO survey (Gilmore et al. 2012). 

We also attempted to add another cluster, NGC 2264, 
to our sample using CoRoT-based rotation periods de- 
termined by Affer et al. (2013). We used these measure- 
ments in an attempt to overcome significant blending in 
'TESS caused by the small angular extent of this clus- 
ter. We combined these periods with vsini data from 
Jackson et al. (2016), but we did not perform further 
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analysis after determining that NGC 2264 would not 
yield enough sini values (2 10-15, Jackson & Jeffries 
2010; Corsaro et al. 2017) to constrain the inclination 
distribution. 


3. ANALYSIS 


To determine projected inclinations and analyze their 
distribution in each new cluster (along with a re-analysis 
of the previously studied clusters for uniformity), we fol- 
lowed the steps described in Papers I and II with minor 
changes. The following subsections list the techniques 
we used to determine each of the components of Equa- 
tion 1. 


3.1. Star Selection 


We began by identifying and removing CG18 members 
whose Gaia G magnitude and Gp — Grp color placed 
them on the equal-mass binary main sequence of their 
cluster's color-magnitude diagram. In addition, we re- 
moved binary and pulsating stars identified by SIMBAD 
(Wenger et al. 2000). Section 3.1 of Paper I describes 
the motivation for removing equal-mass binaries to avoid 
systematic errors in a cluster's inclination determina- 
tions. We avoid pulsating stars in order to minimize 
confusion between pulsation and rotation periods. 

We also removed stars showing excessive astrometric 
noise with Gaia, as quantified by having a Reduced Unit 
Weight Error (RUWE) value greater than the 99th per- 
centile of values estimated to represent single stars (e.g. 
Section 3.1 of Paper II, Belokurov et al. 2020). In order 
to optimize the number of reliable resulting inclination 
determinations for each cluster, we established a CG18 
membership probability threshold at Pmem > 0.5 for 
the final part of the analysis. This threshold balanced 
our preferences for a sample of stars as large as possible 
while containing only probably cluster members. 


3.2. Radius Values 


We used the isochrones software (Morton 2015) to 
fit MIST stellar evolution models (Choi et al. 2016) to 
each star's SED. We built SEDs from the catalogs listed 
in Section 2. We applied an estimated correction to the 
stars’ Gaia EDR3 parallaxes using zero-point (Linde- 
gren et al. 2020). 

We referenced 3D Bayestar19 dust maps (Green et al. 
2019) to estimate reddening for clusters with decli- 
nations greater than —30°. We supplemented these 
maps with a combination of 2D maps from Schlafly & 
Finkbeiner (2011) and targeted studies. We list the full 
set of Gaussian priors we used for SED fitting along with 
references in Table 2. 

In some cases, we assigned the priors a greater uncer- 
tainty than quoted in individual studies in order to avoid 


over-constraining the resulting radius values. Due to the 
lack of a published age estimate for Gulliver 6, we set 
the same prior as ASCC 19 owing to the clusters’ prox- 
imity, their potential to be part of an aggregate group 
of clusters (e.g. Piecka & Paunzen 2021; Soubiran et al. 
2018), and their sharing of one member star according 
to CG18. We eliminated this star from further analysis 
due to its unknown constituency. 

After completing SED fitting, we identified stars 
whose nested sampling global log-evidence was less than 
a threshold we determined empirically for each cluster 
(see Section 4.1 of Paper II). This had the effect of re- 
moving 1-596 of stars in each cluster with the poorest 
SED fits from further analysis. 


3.3. Rotation Periods 


To determine rotation periods, we used a combina- 
tion of TESS light curves from the PATHOS project 
(Nardiello et al. 2019, 2020; Nardiello 2020; Nardiello 
et al. 2021) and eleanor software (Feinstein et al. 2019). 
As in Paper I, we generated each light curve's autocorre- 
lation function (ACF, e.g. McQuillan et al. 2013, 2014) 
and periodgram (e.g. Nielsen et al. 2013) to provide pe- 
riod estimates. Combining the information from these 
analysis tools with a Gaia-based list of nearby blended 
stars and prior expectations of color-period relations in 
open clusters (e.g. Barnes 2003; Kovács 2015; Angus 
et al. 2019), we manually vetted each star in the sam- 
ple. We did not apply a membership probability thresh- 
old for this part of the analysis in order to maximize the 
number of reported periods for stars that are potentially 
cluster members. 

Our evaluations yielded periods for some stars in each 
cluster that demonstrate light curve modulation consis- 
tent with rotating starspots (as opposed to oscillations, 
pulsations, or eclipses). The stars passing the vetting 
process are also not blended with other stars of similar 
color or brightness, which would confuse the assignment 
of a periodic signal to the correct star. Finally, we ap- 
plied extra scrutiny to signals near harmonics of the ~ 
13.7-day TESS satellite's orbital period due to the in- 
fluence of scattered light at those periods. 


3.4. Projected Rotation Velocities 


After accessing the various data sources listed in Ta- 
ble 1, we referenced the source papers to determine an 
appropriate threshold for vsini values. We do not in- 
clude measurements below the threshold in our inclina- 
tion analysis. Section 3.3 of Paper II further discusses 
the motivation of excluding projected rotation velocity 
measurements too close to the predicted macroturbulent 
velocity of a star (see also Kamiaka et al. 2018). 
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Table 1. Clusters in this study (with alternate names), galactic longitude and latitude, Gaia EDR3 distance, mean ages and 
standard deviations from SED fitting posteriors, mass estimates, sources of vsin£, and vsini threshold established to select 


reliable measurements. 


Cluster UL] b [°] d [pc] Age [Myr] M [Mo] vsini src. | Thresh. (km s^!) 

Collinder 69 (A Orionis cl.) | 195.118 | -12.070 392 + 10 4.9 + 0.9 150-250 [1,2] 9 10 
ASCC 16 201.204 | -18.375 | 347.4 + 7.5 (1.3 2:2 ~ 250 [1] 9 10 
ASCC 19 204.884 | -19.371 | 355.5 + 7.8 8.0 + 5.4 ~ 400 [1] 9 10 
Gulliver 6 205.249 | -18.139 415 + 12 8.6 + 5.8 = 9 10 
Pozzo 1 (y Velorum cl.) 262.727 | -7.713 | 346.6 + 5.9 | 10.2 + 3.7 | 400-2000 [3,1 10 10 
BH 56 (vdB-Ha 56) 264.457 1.582 910 + 42 11.8 + 3.8 ~ 150 [1] 9 10 
NGC 2547 264.426 | -8.605 | 386.8 + 6.7 | 31.7 + 7.5 50-450 [1,4] 11 5 
Alpha Per cl. (Melotte 20) | 147.350 | -6.324 | 174.4 + 2.1 | 58.3 + 4.9 | 350-7000 [5,1 9 10 
Blanco 1 14.906 | -79.304 | 236.3 + 4.4 105 + 15 300-3000 [6,1 12 5 
NGC 2422 (M47) 230.961 3.134 476 + 11 155.9 + 8.2 | 200-450 [1,7] 13 5 
NGC 2548 (M48) 227.841 | 15.377 772 + 31 455 + 25 370-1100 [8,1 14 12 

1 Piskunov et al. 2008, ? Bayo et al. 2011, ? Jeffries et al. 2009, 4 Jeffries et al. 2004, 5 Sheikhi et al. 2016, 6 Zhang et al. 2020, 7 Prisinzano 


et al. 2003, $ Balaguer Núñez 2006, ? SDSS DR17 (Abdurro'uf et al. 2022), 10 Jeffries et al. 2014, * Jackson et al. 2016, 12 Mermilliod 
et al. 2009, 13 Bailey et al. 2018, 14 WIYN Open Cluster Study (Sun et al. 2020) 


For SDSS APOGEE-2 data, we established a 10 km 
s71 vsini threshold and assumed uncertainties of ~ 2 
km s^! or 10%, whichever was greater (e.g. Tayar et al. 
2015; Simonian et al. 2020). With the above thresh- 
old, this meant that every APOGEE vsini value we 
analyzed had a 10% uncertainty. To further exclude 
poor-quality data, we removed stars having the fol- 
lowing SDSS quality flags: (ASPCAPFLAGS): VSINI BAD, 
SN BAD, STAR BAD, CHI2 BAD, VMICRO BAD; (STARFLAGS): 
BAD PIXELS, VERY BRIGHT NEIGHBOR, LOW SNR. 

In two other cases, we adopted a vsini threshold 
greater than suggested by source papers. We set thresh- 
olds of 5 km s^! for Blanco 1 data (Mermilliod et al. 
2009) (consistent with our analysis of the same source 
in Paper II) and NGC 2422 (Bailey et al. 2018). 

For Pozzo 1 data (Jeffries et al. 2014), we estimated 
vsini uncertainties by comparing signal-to-noise ratios 
(S/N) for the cluster's spectroscopy with those reported 
by Jackson et al. (2016) for NGC 2264, NGC 2547, and 
NGC 2516. These latter data were collected with the 
same GIRAFFE spectrograph that was used for the ma- 
jority of Pozzo 1 spectroscopy by Jeffries et al. (2014). 
We computed the mean reported uncertainty for Jackson 
et al. (2016) vsin? data within 18 bins of S/N. We then 
performed linear interpolation to estimate the v sini un- 
certainties for Pozzo 1 as a function of S/N using the 
binned data. Above the cluster's reported threshold of 
10 km s~!, we estimated a mean uncertainty of ~ 20% 
in v sin i. 


3.5. Accounting for Systematic Error 


For each cluster, we used the technique described in 
Section 3.3 of Paper II to remove v sin? values suffering 
from potential systematic error. Figure 1 shows the pe- 
riod and effective temperature Tig values we determine 
for stars having vsin? values passing the above criteria 
for selection. For each cluster, we used a MIST isochrone 
generated from the cluster's mean parameters to assign 
radius values to a grid of P and Teg values spanning the 
plotted parameter space. For each point on the grid, we 
calculated an equatorial rotation velocity (i.e. a maxi- 
mum usin?) from the associated P and R values. The 
red contour in each plot represents combinations of pa- 
rameters that result in a maximum vsin i equal to our 
reliability threshold for that cluster. 

Since all stars in Figure 1 already survived the thresh- 
old process, points located above the contour have a 
vsini measurement greater than the limit set by the 
star's rotation period and radius (after allowing for ob- 
served scatter in these two parameters). The color- 
coding of the panels indicates that these stars’ prelimi- 
nary sin i values are often greater than unity, emphasiz- 
ing the potential presence of a systematic error in v sin i. 
We removed stars falling above their cluster's contour to 
decrease the influence of systematic errors on the result- 
ing inclination distributions. In addition, we removed a 
small number of stars having a highly discrepant radial 
velocity compared to the typical values for their cluster 
(marked by a red diamond). 
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Table 2. Gaussian prior values for isochrone fitting and their sources. The three priors are the cluster age in Myr, Fe/H 


metallicity in dex and V-band reddening in mag. 


Cluster Parameter Value Source 
age 5.2+1.2 Bayo et al. (2011); Barrado y Navascués et al. (2004) 
Collinder 69 Fe/H —0.169 + 0.093 Abdurro’uf et al. (2022) 
Av 0.49 + 0.22 Green et al. (2019) 
age 10 + 2.3 Bossini et al. (2019); Kharchenko et al. (2013) 
ASCC 16 Fe/H]  —0.092 + 0.087 Abdurro'uf et al. (2022) 
Av 0.188 + 0.065 Green et al. (2019) 
age 20.0 + 9.5 Bossini et al. (2019); Kharchenko et al. (2013) 
ASCC 19 Fe/H —0.13 + 0.12 Abdurro'uf et al. (2022) 
Av 0.49 + 0.28 Green et al. (2019) 
age 20.0 + 9.5 same prior as ASCC 19 
Gulliver 6 Fe/H —0.15 + 0.11 Abdurro'uf et al. (2022) 
Av 0.64 + 0.36 Green et al. (2019) 
age 20.0 + 4.6 Jeffries et al. (2017, 2009) 
Pozzo 1 Fe/H —0.04 + 0.05 Jeffries et al. (2014); Spina et al. (2014) 
Av 0.131 + 0.055 Jeffries et al. (2009) 
age 17.4 + 4.0 Kharchenko et al. (2005) 
BH 56 Fe/H —0.12 + 0.11 Abdurro’uf et al. (2022) 
Av 0.2 + 0.2 Kharchenko et al. (2005) 
age 34.7 + 4.0 Jeffries & Oliveira (2005) 
NGC 2547 Fe/H —0.03 + 0.06 Magrini et al. (2015) 
Ay 0.37 + 0.16 Naylor & Jeffries (2006) 
age 63.1 + 7.3 Basri & Martin (1999) 
Alpha Per Fe/H —0.08 + 0.14 Abdurro’uf et al. (2022) 
Ay 0.017 + 0.062 Green et al. (2019) 
age 115+11 Gaia Collaboration et al. (2018a) 
Blanco 1 Fe/H —0.09 + 0.13 Buder et al. (2021) 
Ay 0.05 + 0.05 Schlafly & Finkbeiner (2011) 
age 155 + 18 Cummings & Kalirai (2018) 
NGC 2422 Fe/H —0.08 + 0.21 Bailey et al. (2018) 
Av 0.284 + 0.097 Green et al. (2019) 
age 447 + 52 Barnes et al. (2015) 
NGC 2548 Fe/H —0.03 + 0.10 Cui et al. (2012); Hamer (2021) 
Ay 0.145 + 0.062 Green et al. (2019) 


We also estimated the effect of magnetic radius in- 
flation due to spots on the surface of the stars in our 
sample. The fraction of a star’s surface that is covered 
by spots can range from a few to tens of percent, depend- 
ing on age and mass (e.g. Nichols-Fleming & Blackman 
2020; Jackson & Jeffries 2014). Starspots inhibit con- 
vection, reducing a star’s effective temperature. Under 
the simplifying assumption that the star’s luminosity 
can be described by the Stefan-Boltzmann equation, a 
reduction in effective temperature at the same luminos- 
ity will increase the stellar radius (see Sections 5.2 of 
Paper II and 5.6 for further discussions of radius in- 


flation). Since MIST isochrones do not account for the 
presence of starspots, our radius values from SED fitting 
could be underestimated, resulting in an overestimation 
of sini (See Equation 1). This bias could be amplified by 
the requirement of our analysis that stars display light 
curve modulation assumed to be caused by starspots. 
We accessed SPOTS isochrones (Somers et al. 2020), 
which incorporate starspots into their solar-metallicity 
stellar models. The isochrones assume that the effec- 
tive temperature of starspots are 80% that of the non- 
spotted stellar surface. To obtain an estimate of ra- 
dius inflation’s effect without overcompensating, we se- 
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lected isochrones corresponding to the lowest nonzero 
fractional spot coverage modeled, fspot = 0.17. For each 
cluster, we compared MIST-derived radii at the same 
age and solar metallicity interpolated onto the mass val- 
ues of SPOTS. We computed the mean fractional dif- 
ference between the two sets of radius estimates, find- 
ing that the radii we report may be underestimated by 
~ 2 — 1396. To incorporate this prediction into our 
analysis without excessively modifying our results, we 
treated the mean fractional difference in radius as a sys- 
tematic error and added it in quadrature to the lower 
error bar of each sin i value. We only did this for stars up 
to 1.3Mo to avoid extrapolating the SPOTS isochrones 
above their available range of masses. 

Finally, following Section 5.6 of Paper II, we estimated 
the effect of differential rotation on sin ? values, which we 
summarize here. Starspots at nonzero latitudes on the 
stellar surface will rotate more slowly that spots at the 
equator, leading to overestimated rotation periods (e.g. 
Hirano et al. 2014). In Paper IL, we predicted that if 
starspots exist at solar-like latitudes (1 = 20+ 20?) with 
a solar-like rate of differential rotation (a fractional dif- 
ference between equatorial and polar rotation rates be- 
tween 0.1 and 0.2), the resulting rotation periods will be 
overestimated by an amount between ~ 1-296 on aver- 
age. 

Projected rotation velocities will be underestimated 
in the presence of differential rotation, since some parts 
of the star are rotating more slowly than the equatorial 
rotation rate. This effect will manifest on the scale of 
~ 5-10% in solar-like rotators (Hirano et al. 2014). Av- 
eraging the combination of opposing systematic effects 
due to differential rotation on rotation periods and v sin 2 
values, we compute a mean systematic underestimate of 
sini by 5.796. We added this fractional error in quadra- 
ture to the upper error bar of our reported inclination 
values. 


3.6. Inclination Distributions 


We combined R, P, and vsin i values to determine the 
sine of each star's inclination and its uncertainty via the 
Bayesian method described by Masuda & Winn (2020). 
We assumed that each parameter was represented by a 
Gaussian uncertainty distribution except for the vsin? 
values of NGC 2547. For this cluster, we instead mod- 
eled each projected rotation velocity as a t-distribution 
with two degrees of freedom to better describe the ex- 
tended tails of each measurement's distribution accord- 
ing to Jackson et al. (2015). 

We modeled the resulting inclination distributions us- 
ing the cone model created by Jackson & Jeffries (2010). 
'This framework imagines stellar spin-axes uniformly dis- 


tributed within a cone having a spread half-angle \ and 
mean inclination a. Isotropic spin-axis distributions 
(A ~ 90°) reduce the distribution to one that is uniform 
in cost, while tightly aligned distributions will have A 
values closer to 0°. Mean inclinations range from pole- 
on orientations (a ~ 0°) to edge-on (o ~ 90°). We also 
employ our modification of this framework to include a 
third parameter representing a fraction of stars aligned 
amid an otherwise isotropic distribution. More details 
about this three-parameter model are offered in Section 
3.5 of Paper II. 

We forward-modeled the effect of establishing a 
threshold in v sin ? on the inclination distributions. This 
selection effect will generally eliminate stars that are ei- 
ther intrinsically rotating slowly or have low projected 
rotation velocity due to a nearly pole-on inclination. 
'This second category of stars is especially important to 
consider, since the preferential elimination of data as an 
implicit function of the quantity of interest will directly 
affect our overall results. As described in Section 3.4 
of Paper II, we used a Monte Carlo simulation to cre- 
ate a population of ~ 10° stars matching the color dis- 
tribution for each cluster's stars. We used a Gaussian 
Mixture Model clustering algorithm to identify rotation 
sequences in the color-period space of our data, and we 
designated a proportionate fraction of simulated stars to 
each sequence based on the results of the clustering. 

We fit polynomial splines to each identified color- 
period sequence to assign rotation periods to simulated 
stars. We used a MIST isochrone to assign radius values 
as a function of color. Combining the simulated radius 
and period values, we computed rotation velocities for 
our simulated population. We paired each rotation ve- 
locity with a sin? value generated from the Jackson & 
Jeffries (2010) cone model with the desired A and a val- 
ues. Thus, we created a collection of model v sin į values 
for which we set the threshold chosen for each cluster. 
By inducing a change in the associated sin i distribution 
after discarding simulated stars below the v sin i thresh- 
old, we modeled the expected selection effect on our data 
after setting the same threshold. 

To compare our sin? determinations for each clus- 
ter with the simulated values generated using the cone 
model, we computed cumulative distribution functions 
(CDFs) for both. We modified the resulting model 
CDFs to account for random uncertainties in R, P and 
vsini. As in Paper II, we assumed normal distributions 
of fractional errors in P and R, and log-normal distri- 
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Figure 1. Period-effective temperature plots for the clusters newly analyzed in this paper. Stars are color-coded by the sin i 
value their parameters yield. Red curves show the contour on which a star having a threshold vsin? value yields sin? = 1, 
incorporating scatter in the P and R values. Above this contour, stars have nonphysical inclinations and are discarded for their 
potential contribution of systematic error to a cluster's inclination distribution. Red diamonds indicate stars removed from the 
sample due to having an outlying radial velocity value indicative of a potential binary system. 
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butions of fractional errors in v sini (e.g. Kovacs 2018).! 
Setting a vsin1 threshold increases the sin? value asso- 
ciated with a given cumulative probability in the model 
CDFs, shifting them rightward compared to a model 
without a threshold. At any given threshold level, the 
magnitude of the shift is proportional to a sin? value’s 
proximity to 0°. Using these models, we performed a 
10,000-step MCMC analysis of parameter space using 
the same Gaussian likelihood function as Equation 2 of 
Paper II: 


InLF = —0.5 * > E (2102) 4- pc . (2) 


i=l 


Here, n is the number of sini values for a cluster, x; 
is a determination of sini, v; is the uncertainty, and £; 
is the model's predicted inclination at the same cumula- 
tive probability as z;. The number of degrees of freedom, 
v =n — m, is used to calculate the reduced chi-squared 
value x% a, having either m = 2 or m = 3 fitted param- 
eters (a, A, and f for the three-parameter model). For 
the two-parameter model, as in Paper II, we discarded 
a 1000-step burn-in before marginalizing the posteriors 
to obtain probability distributions (PPDs) for o and A. 
To facilitate the use of asymmetric error bars, we used 
either the upper or lower value of c; depending on the 
sign of the difference between data and model sin? val- 
ues. For x; — £; > 0, we used the lower value of c;, and 
vice-versa. We also explored a three-parameter model 
when two parameters did not provide a satisfactory fit 
to the data (see Section 4.3). 


3.7. Cluster Kinematics 


'To search for connections between stellar rotation and 
bulk motion within a cluster, we again computed plane- 
of-sky and line-of-sight (LOS) internal kinematics based 
on Kamann et al. (2019) and Paper I (Section 3.7) and 
Paper II (Section 3.6). 

To calculate plane-of-sky kinematics in multiple bins 
of radial distance from cluster center, we accounted for 
bias in the Gaia EDR3 proper motions of bright stars 
(Cantat-Gaudin & Brandt 2021), transformed to Carte- 
sian coordinates (Gaia Collaboration et al. 2018b), sub- 
tracted contributions due to each cluster's motion along 
the LOS (van de Ven et al. 2006), and transformed again 
to polar coordinates (van Leeuwen et al. 2000). In this 


1 While we assumed that each individual v sin íi measurement was 
quantified by either a normal or t-distribution depending on the 
measurement source, we found that the overall ensemble of frac- 
tional vsin? uncertainties for each cluster was better described 
by a log-normal distribution. 


coordinate system, a star's radial unit vector points in 
the direction opposite that of cluster center, while the 
tangential unit vector is oriented perpendicular to the 
radial direction and is positive in the clockwise direc- 
tion on the sky (see Footnote 5 of Paper IT). 

To analyze LOS kinematics, we used Gaia EDR3 LOS 
velocities for Blanco 1 and NGC 2422. For all other 
clusters lacking a sufficient quantity of Gaia LOS ve- 
locities, we consulted the targeted spectroscopic studies 
listed in Table 1. After discarding likely spectroscopic 
binary outliers from the primary LOS velocity distri- 
bution, we determined mean LOS kinematics for each 
cluster. The smaller number of LOS measurements com- 
pared to proper motions prevented us from drawing sig- 
nificant insights from binning the results by radial dis- 
tance. 


4. RESULTS 


We plot color-magnitude diagrams for each cluster in 
Figure 2. The black points in these figures represent 
all stars on the CG18 membership list. In each panel, 
we also plot the subsets of cluster members with rota- 
tion period measurements and sini determinations. The 
following subsections and their associated tables and fig- 
ures further describe our results for the 11 clusters we 
analyzed. 


4.1. Rotation Periods 


For the 11 clusters, we report 1618 rotation periods 
we determined with TESS in Table 3. We plot rota- 
tion periods versus Gpp — Grp color for all clusters in 
Figure 3. For clusters older than 30 Myr, we include 
an empirical isochrone for the slow-rotating sequence of 
stars (red dashed line) based on either Barnes (2003) 
or Angus et al. (2019). For the former, we estimated 
Gpp — Grp colors as a function of B — V using Gaia 
documentation.” The gyro isochrones serve to aid the 
visualization of each cluster’s rotation sequence and sup- 
port the ages and reddening values that result from our 
SED fitting analysis. We determined rotation periods 
for between ~ 24-58% of stars across the entire sample, 
with a mean of 33%. 


4.2. Inclinations 


We determined sini values for a total of 216 stars 
across the 11-cluster sample. Table 4 lists these values 
along with their upper and lower uncertainties, includ- 
ing the systematic errors described in Section 3.5. The 
mean fractional uncertainty of the sin i determinations 


? *Photometric relationships between Gaia photometry and exist- 
ing photometric systems," C. Jordi 
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is ~ 1396. Table 4 also lists values for Tog and R from 
SED fitting, P from our rotation period analysis, and 
usin? from the corresponding data source (see Table 1). 

Figure 4 shows histograms for each cluster's inclina- 
tions compared to an isotropic distribution. Figure 5 
shows plots of the cone model posteriors (normalized by 
the maximum value) in « and A parameter space along 
with 1D marginalized PPDs. Figure 6 plots the empiri- 
cal CDFs for each cluster's inclinations, an isotropic dis- 
tribution, a CDF corresponding to the PPD-peak values 
of A and a, and CDFs representing model parameters 
drawn from the 6896 confidence intervals of the PPDs. 
Table 5 summarizes the PPD-peak and 9596 confidence 
intervals in A and o for each cluster's two-parameter 
MCMC run, along with the reduced x? and maximum 
natural logarithm of the posterior values for the most 
probable results. 

The numerical simulations of Corsaro et al. (2017) 
suggest that stars below 0.7 Mo do not show spin align- 
ment even if more massive cluster members do. We ran 
a second iteration of our MCMC analysis on the stars 
greater than this mass threshold for each cluster, find- 
ing no significant differences in results compared to the 
entire-sample analyses. 


4.3. Three-parameter Model Runs 


We performed a grid search of 3-parameter model 
space for each cluster in order to identify significant 
improvements over two-parameter model fits. The f 
parameter of this model quantifies the fraction of stars 
aligned in a cone described by a mean inclination Qaniso 
and spread half-angle Aaniso. The remaining 1 — f frac- 
tion of stellar spins are modeled as isotropic. We com- 
puted the natural log-posteriors for values of a and A 
between 10? and 90? (spaced in 10? intervals) and for f 
between 0.1 and 0.9 (spaced by 0.1). 

We plot two visualizations of this grid search in Fig- 
ures 7 and 8. Each panel of these figures represents a 
chosen value of f. The posterior plots in each panel 
highlight the most probable values of Qaniso and Aaniso- 
The more similar each posterior plot looks to the two- 
parameter version shown in Figure 5, the less distinctive 
that three-parameter fit is from the two-parameter fit. 

In one case, we noticed a significantly distinct three- 
parameter model solution and ran the full MCMC anal- 
ysis to determine the optimal fit. We further exam- 
ined the inclinations of stars in NGC 2548, since the 
grid search identified better three-parameter fits than 
the best-fitting two-parameter model, which yielded the 
greatest reduced x? value (6.6) among the clusters an- 
alyzed. Allowing for a subset of aligned stars, we de- 
termined a maximum-posterior three-parameter model 


solution of Qaniso = 24°, Aaniso = 4°, and f = 0.3. The 
reduced xy? value corresponding to this parameter com- 
bination is 3.4, and the log-posterior is 33.3. Figure 9 
compares this cluster's two- and three-parameter model 
fits. 

For NGC 2548, we selected six stars of the 22 cor- 
responding to a mean inclination and fraction of stars 
consistent with our three-parameter model expectation. 
The inclinations of this group modeled as aligned are 
marked as blue in Figure 9, while we label the remaining 
stars as an isotropic subset. As in Paper II, we studied 
the positions and proper motions of the selected stars to 
search for groupings in their position angles and plane- 
ofsky positions (Figure 10) as well as internal proper 
motions (Figure 11). 

A Kuiper test on the position angles of the se- 
lected group under the null hypothesis of random an- 
gles yielded a p-value of 0.08. A bootstrapped, two- 
sample, two-dimensional K-S test (Peacock 1983; Fasano 
& Franceschini 1987; Press et al. 2007)? comparing 
plane-of-sky positions for the selected group to those 
of the background group found a p-value of 0.09. An- 
other two-sample 2D K-S test, this time comparing the 
internal proper motions of the selected and background 
groups, returned a p-value of 0.49. 


4.4. Internal Kinematics 


Figure 12 plots our results for the binned internal kine- 
matics analysis of each cluster, and Table 6 lists our de- 
terminations of mean kinematic parameters in the tan- 
gential, radial and LOS directions. As in Paper II, we 
have subtracted the predicted trend due to motion along 
the LOS from the radial values. 


5. DISCUSSION 


In this section, we begin by discussing our interpreta- 
tions of the results for the individual clusters in separate 
subsections. Since a small sample of inclinations taken 
from a random distribution can appear to be anisotropic 
(as quantified in simulations by Corsaro et al. (2017)), 
we only included clusters having 15 or more sin i values 
in our ensemble analysis. As a result, we do not dis- 
cuss the inclination distributions of ASCC 16, Gulliver 
6, Pozzo 1, and BH 56 in this section. 

We continue with a discussion of the added inferences 
we make from the full constrained sample among the 15 
clusters in this work and Papers I and II. We obtain en- 
semble insights using two different approaches: classifi- 
cations of spin-axis distributions on a cluster-to-cluster 


3 ndtest Python code written by Zhaozhou Li, available at https: 


/ / github.com/syrte/ndtest 
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Table 4. Relevant quantities and determinations of the inclinations of open cluster members. 
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NoTE—Table 4 is published in its entirety in the machine-readable format. This table includes updated data for NGC 2516, 


Pleiades, Praesepe, and M35 data (Papers I and II) due to the revised procedures detailed in this paper. 
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Figure 2. Color-magnitude diagram plots for all clusters. Orange triangles indicate stars for which we measured rotation 


periods with TESS. Cyan diamonds show stars with sini values determined by this work. Yellow star symbols indicate reddened 
and extincted solar values at each cluster’s distance. 
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Table 5. Inclination distribution parameters based on the Jackson & Jeffries (2010) two-parameter cone model. We list the 
number of rotation period and sin? values for each cluster, the values corresponding to PPD peaks for the spread half-angle A 
and mean inclination a, and their 95% confidence ranges. We also report the reduced x? and logarithm of the unnormalized 
posteriors corresponding to the PPD-peak parameters for each cluster. The final column lists the classification corresponding to 
the best-fitting model for each cluster’s spin-axis distribution: isotropic/moderately aligned (I), unconstrained (U), or partially 


aligned (P). 


Cluster Npe | Nan; | A PPD peak (95% conf. range) | a PPD peak (95% conf. range) | x24, | In(post.) | Class. 
Collinder 69 | 202 26 86° (A > 22? 54° (5? <a < 62°) 1.1 24.2 I 
ASCC 16 | 129 | 10 86° ( > 6°) 88° (16° < a < 88°) 1.0 9.1 U 
ASCC 19 157 17 89° (A > 15° 53° (9° <a < 85°) 0.6 18.0 I 
Gulliver 6 | 121 | 10 59° (25° < A < 81°) 37% (4° < a < 52°) 08 | 123 U 
Pozzo 1 110 13 90° (A > 11? 67? (13° < a < 87°) 0.6 11.4 U 
BH 56 45 | 11 89° (A > 5°) 85° (19° < a < 88°) 0.9 9.8 U 
NGC 2547 79 27 89° (A > 45° 51? (9° <a < 79°) 1.7 45.1 I 
Alpha Per 294 30 89° (A > 6°) 69° (12° < a < 87°) 0.4 25.0 I 
Blanco 1 187 29 86° (A > 21° 57° (9? <a < 85°) 1.3 31.5 I 
NGC 2422 159 21 50° (11° <A « 71?) 73° (62° <a < 88°) 0.5 39.2 U 
NGC 2548 135 22 58° (35° < A < 74°) 31° (10° < a < 45°) 6.6 30.6 P 
NGC 2516 | 158 | 44 90° (A > 34°) 53° (9° <a < 84°) 0.2 | 72.0 I 
Pleiades 159 39 90? (A > 29?) 53° (9° <a < 83°) 0.4 55.6 I 
Praesepe 947 36 87° (A > 28°) 52° (9? <a < 83?) 0.8 41.8 I 
M35 1146 | 69 59° (37° < À < 81°) 36° (6° <a < 50°) 2.3 | 100.6 P 
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Table 6. Mean determinations of internal kinematics for each cluster. 


Cluster Direction 6. [°] vo [km s^ !] vrot [kms !] | e [kms !] 
Tangential — —0.03' + 0.033 — 0.685 + 0.024 
Collinder 69 Radial m 0.842 + 0.046 — 0.953 + 0.033 
LOS 110+10  25.09+0.14  0.90+0.23 1.162+ 0.086 
Tangential — —0.115 + 0.028 — 0.452 + 0.021 
ASCC 16 Radial — 0.278 + 0.028 — 0.445 + 0.020 
LOS 2002-40 20.33 +0.24 0.29+9:23 0.83+5:09 
Tangential - —0.065 + 0.030 mE 0.502 + 0.022 
ASCC 19 Radial — 0.485 + 0.029 — 0.500 + 0.021 
LOS 340440 19.88+0.26  0.67+0.35 — 1.544 0.19 
Tangential = —0.079 + 0.052 = 0.739 + 0.038 
Gulliver 6 Radial — 0.268 + 0.051 — 0.762 + 0.037 
LOS 260490 32.47 + 0.21 0.251927 1.021014 
Tangential — —0.021 + 0.027 — 0.496 + 0.019 
Pozzo 1 Radial — 0.042 + 0.028 — 0.526 + 0.020 
LOS 1102-90 17.811+0.085 — 0.09070 004 — 0.651 + 0.065 
Tangential — —0.22 + 0.13 — 0.89 + 0.10 
BH 56 Radial — 0.25 + 0.15 — 1.01 +0.11 
LOS 1004-30 23.49+0.31 0.864048 1.29 + 0.24 
Tangential — —0.080 + 0.033 — 0.431 + 0.024 
NGC 2547 Radial — 0.017 + 0.034 — 0.446 + 0.025 
LOS 3102-50 11.820+0.086 0.16+0.12 — 0.679 4- 0.065 
Tangential — —0.491 + 0.038 = 0.902 + 0.027 
Alpha Per Radial — —0.290 + 0.041 = 0.958 + 0.029 
LOS 070 —1.814 0.073 — 0.090*0093 0.704 + 0.056 
Tangential — 0.064 + 0.030 — 0.472 + 0.021 
Blanco 1 Radial — 0.112 + 0.030 == 0.485 + 0.022 
LOS 180110 — 4.31 0.24 0.241926 Witt ae 
Tangential — 0.016 + 0.028 — 0.531 + 0.021 
NGC 2422 Radial — 0.001 + 0.029 — 0.553 + 0.021 
LOS 10+70 36.00 0.58 use C632 50 
Tangential — 0.066 + 0.036 — 0.662 + 0.026 
NGC 2548 Radial — 0.078 + 0.038 — 0.687 + 0.027 
LOS 090 9.75 + 0.18 Ü 2d 0.90 + 0.16 


NOTE—4. is the position angle of each cluster's LOS rotation axis; vo is the velocity in the specified direction, Vrot is the LOS 
rotation velocity, and o is the velocity dispersion in the specified direction. 


scale, and allocations on a star-by-star scale based on 
model fitting results. 


5.1. Newly Analyzed Clusters 
5.1.1. Collinder 69 


This cluster’s CDF is shifted rightward (towards 
greater values of sin ?) compared with all possible model 
CDFs. Despite this shift, isotropic spins are preferred by 
the spread half-angle PPD. The Collinder 69 inclination 
distribution is likely influenced by magnetic radius infla- 
tion due to the young age of the cluster. The systematic 


shift may also be caused an overestimate of SDSS DR17 
projected rotation velocities or TESS rotation periods. 

Despite the visible offset between the model and data 
CDFs, we can still make statistical inferences about the 
alignment or isotropy of the cluster’s spin-axis distribu- 
tion, since we have made informed efforts to account for 
inflated radii and uncertainties in APOGEE vsini val- 
ues, as well as the influence of differential rotation on 
period measurements. By increasing (or establishing) 
the appropriate uncertainties on the data, we decrease 
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Figure 3. Color-rotation period plots for all clusters. All plotted periods are determined with TESS data by this work. For 
clusters older than 30 Myr, we plot an empirical sequence of slow rotators (Barnes 2003; Angus et al. 2019) which validate the 


parameters we determined via SED fitting. 


the statistical significance of systematic errors at the ex- 
pense of precision in the results. 

Nonetheless, we find that the spread half-angle for 
Collinder 69 has a meaningful constraint, with A > 22° 
disfavoring tight alignment at 9596 confidence. The pos- 
terior plot for this cluster displays the familiar *penin- 
sula" shape we first described in Section 4.2 of Paper II: 
isotropic and near-isotropic regions of parameter space 


(A Z 65?) are favored, extending to moderate align- 
ment at moderate mean inclinations (30° < A < 65°; 
35° S a S 75°). These regions of the posterior plot 
are surrounded by lower-probability parameter combi- 
nations. Similar to our results for NGC 2516 (Paper 
I) and the Pleiades and Praesepe (Paper II), we cannot 
distinguish between isotropy and moderate alignment 
for Collinder 69. However, we make further inferences 
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Figure 4. Histograms of projected inclinations for all clusters (blue), with isotropic distributions overlaid in orange. Error bars 


assume Poisson statistics. 


about this and other clusters showing similar peninsula- 
shaped posterior plots in Sections 5.3 and 5.4. 


5.1.2. ASCC 19 


This cluster’s posterior plot displays a peninsula shape 
that is comparable to that of Collinder 69. The log- 
posterior of the PPD-peak model is smaller for ASCC 
19, reflecting the smaller number of sin values for this 
cluster. The reduced x? value for the PPD-peak model 


is also smaller than Collinder 69, reflecting the better 
CDF fit for ASCC 19 and less systematic offset among 
inclination values. 

The optimal model’s value of xX2,, < 1 indicates the 
possibility of overfitting or an overestimate of uncertain- 
ties. Because the 2D cone model has low complexity in 
using two parameters to model the CDF, overestimated 
uncertainties may be the more likely source of the re- 
duced y? value. Since vsin i uncertainties are the dom- 
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Figure 5. Two-parameter cone model log-posterior plots for all clusters, normalized by each cluster’s maximum log-posterior 
value. Each panel also displays marginalized PPDs for a and A. Each plot shows a randomly selected 2% of the full MCMC 


chains. 


inant source of error in each sin determination and we 
used an approximate estimate for the APOGEE vsin 
uncertainties (10% for all values above 10 km s^), it 
is possible that we did not accurately represent the un- 
certainty for every value. Absent further analysis of the 
APOGEE spectra, our uncertainty estimates still lead 
to informative constraints on spin alignment with a de- 
termination of A > 15? at 95% confidence. Overall, 


we classify ASCC 19 as either isotropic or moderately 
aligned. 


5.1.8. NGC 2547 


The PPD-peak result for NGC 2547 yields a greater 
log-posterior than for any of the previous clusters. The 
posterior plot forms a peninsula shape similar to those 
of Collinder 69, ASCC 19, and Pozzo 1, but in this clus- 
ter's case the extent of high-probability parameter space 
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Figure 6. Inclination CDF plots for all clusters. Solid red curves correspond to PPD peaks for a and A. Dashed curves show 
the isotropic distribution. Pink CDFs are generated with cone model parameters drawn from the 68% confidence interval for 
each cluster. 
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Figure 7. Blanco 1 log-posterior plots for a grid of three-parameter model values. Values of Ganiso and Aaniso range from 10° to 
90° in 10° increments, and f ranges from 0.1-0.9 in increments of 0.1. Each panel represents the results for specific value of f. 
For low values of f, the aligned subset’s parameters are not constrained, and small amounts of systematic error can be overfit 
by the more complex model. As f increases, the posterior plots approach their two-parameter counterpart, degenerate between 
isotropic and moderately aligned (see Figure 5). 
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Figure 8. Same as Figure 7, but for a cluster not well-described by isotropy or moderate alignment (NGC 2548). There is 
a significant preference for a specific combination of three-parameter values. No panels are comparable to the two-parameter 
posterior plot for this cluster (see Figure 5) because the corresponding posteriors are much smaller than the optimal three- 


parameter model fit. 


is smaller. The reduced x? value for NGC 2547 sug- 
gests a potential underfit or underestimated uncertain- 
ties, which may be the source of the tighter constraints. 
Since we still see the degenerate solution shared by the 
aforementioned clusters, we add NGC 2547 to the set of 
clusters classified as isotropic/moderately aligned. 


5.1.4. Alpha Per 


This cluster shows a rightward offset of nearly all of 
its sin? values compared with models. As a result, there 
is not a strong constraint on the cluster's spin-axis dis- 
tribution. We have classified other such clusters in our 
sample as having an unconstrained spin-axis distribu- 
tion when the model parameters are loosely constrained. 
However, a previous study of Alpha Per by Jackson & 
Jeffries (2010) determined that the cluster has either 


isotropic or moderately aligned spins (A > 40° at 90% 
confidence and an optimal value of A = 90°). Since this 
scenario is consistent with our result, we adopt it as our 
qualitative interpretation of Alpha Per. Thus, we treat 
the cluster as having cone model parameters degener- 
ate between isotropy and moderate alignment for the 
ensemble analysis. 


5.1.5. Blanco 1 


We classify this cluster as isotropic/moderately 
aligned, since its peninsula-shaped posterior plot 
matches the pattern we have previously identified. The 
three-parameter multi-panel plot in Figure 7 is represen- 
tative not only of Blanco 1 but also the other clusters 
classified as isotropic/moderately aligned. Studying this 
plot in more detail helps to explain differences between 
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Figure 9. Three-parameter cone model fit to NGC 2548 
inclinations (blue) compared to the optimal two-parameter 
model (red, dashed). Blue and black sin i points denote stars 
modeled as aligned and isotropic, respectively. 
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Figure 10. Cluster member positions for NGC 2548 stars 
(gray). Stars with sin į determinations are split into the same 
aligned (blue diamond) and isotropic (orange circle) subsets 
as shown in Figure 9. 


the two- and three-parameter models for the most com- 
monly seen inclination distributions in our sample. 

For low values of f, for which the model includes a 
small fraction of aligned stars, the Qaniso-Aaniso poste- 
rior plot contains approximately uniform high posteri- 
ors for all parameter combinations. This means that 
no constraints can be placed on the orientation, spread 
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Figure 11. Internal proper motions for NGC 2548 stars 
on the single-star main sequence (gray). Stars with sini de- 
terminations are plotted with the same split as described in 
Figure 10. 


half-angle, or even the existence of an aligned subset at 
low values of f. Instead, the three-parameter model uses 
its added complexity to achieve similar or slightly better 
fits than the two-parameter model. We do not interpret 
these results as indicative of an aligned subset of stars 
in these clusters, since the model will tend to overfit to 
slight systematic errors. 

As f increases, the Qaniso-Aaniso posterior plots begin 
to resemble the peninsula shape of the two-parameter 
model MCMC. The background inclination distribu- 
tion of the three-parameter model is always isotropic. 
Thus, if the “aligned” fraction of stars also appears 
isotropic, the resemblance suggests that the three- 
parameter model cannot find a better solution hav- 
ing an aligned subset of stars compared to a single, 
isotropic/moderately aligned inclination distribution for 
all cluster members. 


5.1.6. NGC 2422 


'This cluster displays a signature of anisotropic spins, 
with its PPD-peak model corresponding to moderate 
alignment at a nearly edge-on mean inclination. NGC 
2422 v sin i values (reported by Bailey et al. (2018)) have 
the smallest average fractional uncertainty among all the 
clusters in the sample at ~ 3.3%. These small uncertain- 
ties may either reinforce the preference of alignment over 
isotropy or overconstrain the model to an inaccurate set 
of parameters. The reduced x? value for this cluster is 
suggestive of possible overfitting. 
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Figure 12. Tangential, radial and LOS internal kinematics for all clusters. 
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To further investigate this cluster’s inclinations, we 
compared vsini values reported by Bailey et al. (2018, 
B18) and Jackson et al. (2016, J16) for NGC 2516, a 
cluster we previously analyzed in Paper I. Among the 
14 stars in common between the two studies, the B18 
usin? values are on average 12 + 4% (standard error) 
greater than those of J16. To see how the posteriors 
change for NGC 2422 if this discrepancy is caused by 
a systematic error in B18 vsin? values, we performed a 
new MCMC analysis after adding a 1296 fractional er- 
ror in quadrature to the lower uncertainties of each sin 4 
value. While we do not assume that the J16 values are 
ground-truth by default, we found that if the systematic 
error manifests as an underestimate in J16 v sin? values, 
then the inclinations for NGC 2516 are fit more poorly 
by all possible two-parameter combinations. If we in- 
stead treat the systematic error as an overestimate in 
B18 values, we achieve a better fit for the optimal NGC 
2422 model (as quantified by the reduced x? value). We 
found that our new analysis of NGC 2422 resulted in an 
isotropic PPD-peak value of A = 90°, with the spread 
half-angle constrained to be > 13° at 95% confidence. 

We also experimented with subtracting 1296 from each 
sin? value while keeping the uncertainties unchanged. 
'The resulting MCMC run yielded a similar result as the 
previous one, finding A = 89° and A > 14? at PPD-peak 
and 9596 confidence, respectively. These two additional 
analyses show that while the results we report for NGC 
2422 favor spin alignment, the known discrepancy in 
vsini values reported by B18 and J16 for NGC 2516 is 
enough to change the NGC 2422 result to one that favors 
isotropy. Thus, while we do not change the statistics we 
report for this cluster, we emphasize the potential for 
systematic v sin į error to be the cause of the anisotropic 
signature, and due to this uncertainty we do not include 
NGC 2422 in the ensemble analysis. 


5.1.7. NGC 2548 


'The three-parameter multi-panel plot for NGC 2548 
(Figure 8) represents the appearance of a three- 
parameter combination that achieves a significantly bet- 
ter solution than the two-parameter model (See Figure 
9). As opposed to Figure 7, the familiar peninsula shape 
does not appear for any value of f, nor does the two- 
parameter posterior plot become more visible as f in- 
creases. Instead, the three-parameter model reports low 
log-posterior values for all regions of parameter space 
except for a relatively narrow range of values. We saw 
a similar plot for M35 (Paper II), indicating the insuf- 
ficiency of the two-parameter model to describe these 
clusters' inclination distributions. 


'The two-parameter model MCMC run for NGC 2548 
disfavors isotropy, but the inclination distribution is the 
poorest-fit of our 15-cluster sample with a reduced x? 
value of 6.6. The three-parameter results offer an im- 
provement in modeling values of sin? less than 0.5 but 
still fail to capture the distribution of larger values closer 
to and exceeding unity. Even when using the more com- 
plex model, the lowest achieved x?,¿ value of 3.4 quan- 
tifies the high degree of underfitting compared to other 
clusters in the sample. The statistical tests on the po- 
sition angles, locations and proper motions of poten- 
tially aligned stars did not have significant incompati- 
bility with the null hypothesis of randomness. 

It is worth noting that the v sin? values for both this 
cluster and M35, another cluster for which we found an 
aligned subset as a potential explanation, come from the 
same survey, the WIYN Open Cluster Study. M35 and 
NGC 2548 also saw the two largest fractions of vsin? 
measurements eliminated due to suspected systematic 
error (see Section 3.5 and Figure 1). While this elim- 
ination process removed stars likely to yield nonphys- 
ical sin? values significantly greater than unity, it did 
not target vsinz values that result in low inclinations, 
since it is difficult to disentangle truly pole-on or slowly 
rotating stars from systematic underestimates that are 
still physically feasible. The common origin of the v sin? 
data for this cluster and M35 raises the possibility that 
a shared systematic error may be the source of both 
clusters’ disagreement with the two-parameter models. 

'To unambiguously interpret the results of the three- 
parameter model, the sample size of any subset of cluster 
stars must be great enough to rule out apparent align- 
ment, even if the full set of sin? values in the cluster 
is suitable for constraints on the inclination distribu- 
tion using the two-parameter model. Since only 6 stars 
make up the subset of NGC 2548 inclinations modeled 
as aligned, the same limitations we faced with using the 
two-parameter model to describe ASCC 16, Gulliver 6, 
Pozzo 1, and BH 56 are applicable to NGC 2548 and the 
three-parameter model. While the visibly bimodal sin 2 
distribution of this cluster is noteworthy and merits a 
different qualitative label than the aforementioned clus- 
ters in Table 5, we cannot interpret the results of the 
three-parameter model fit as indicative of spin align- 
ment. Given our model's poor goodness-of-fit, the lack 
of a significant position angle, location, or proper mo- 
tion grouping, and the qualities of the spectroscopic data 
described above, it is plausible that systematic error in 
v sin? caused the apparent alignment. 


5.2. Previously Analyzed Clusters 
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NGC 2516 (Paper I), the Pleiades and Praesepe (Pa- 
per II) all show minimal changes from their previously 
determined inclination distributions. 

For NGC 2516, there is a 33% increase in the number 
of reported sini values because we did not apply an 0.7 
Mo threshold for this version of the analysis (see Sec- 
tion 4.2). These three clusters all have peninsula-shaped 
posterior plots with PPD-peak mean inclination values 
within 1? of each other. Optimal spread half-angle val- 
ues are within 3? of each other. 

M35 saw a slightly larger change from its two- 
parameter model constraints in Paper II. The 9596 con- 
fidence interval in the spread half-angle now extends to 
A = 81? compared to its original value of 67°. This 
change is slightly more favorable of isotropy, and may 
be caused by the standardization of stellar magnitudes 
for SED fitting in this paper. We did not include BV RI 
magnitudes from Thompson et al. (2014) in this work, 
perhaps inducing a slight change in the resulting sin? 
values or their uncertainties. The maximum-posterior 
three-parameter model solution for M35 remains simi- 
lar to the optimal parameters in Paper II, favoring a 
small fraction of stars with tightly aligned spins along 
the LOS. The parameters we determined in this work's 
new analysis are Qaniso = 6°, Aaniso = 9°, and f = 0.06 
(corresponding to ~ 4 stars modeled as aligned, again 
too small a subset to conclude physical alignment). 


5.3. Combined Insights: Cluster-by- Cluster 


Our overall 15-cluster sample of spin-axis distribu- 
tions, including the 11 studied in this work, three in 
Paper II and one in Paper I, establishes patterns and 
permits us to draw ensemble-level insights not obtain- 
able from any single cluster. 


5.3.1. Unconstrained Clusters 


The small number of stars with sini values for ASCC 
16, Gulliver 6, Pozzo 1, and BH 56 precludes us from de- 
termining their spin-axis distributions. In addition, four 
of the five clusters whose spin-axis distributions we clas- 
sify as unconstrained are among the youngest members 
of our sample. As we describe in Section 3.5, these young 
clusters are especially likely to experience radius infla- 
tion due to magnetic activity. The probable presence of 
radius inflation in these clusters, manifesting as an un- 
known bias in sini, further complicates a constraint on 
spin alignment or isotropy for these clusters. 

The fifth unconstrained cluster, NGC 2422, likely ex- 
periences less radius inflation than the other four, but 
the potential systematic v sini error we identify in Sec- 
tion 5.1.6 creates enough uncertainty to similarly pre- 
vent constraints of the cluster’s spin-axis distribution. 


The limitations illustrated by these five clusters empha- 
size that future work on this topic would benefit from 
observation campaigns dedicated to using a multi-object 
spectrograph to yield an increased number of v sini val- 
ues for as many Gaia-based cluster members as possible. 
Combined with the refinement and use of radius inflation 
predictions to mitigate systematic error in this param- 
eter, future studies can provide greater constraints on 
younger and less populous open clusters. 


5.9.2. Clusters Modeled as Partially Aligned 


two clusters benefit from improved fits with the three- 
parameter model compared to insufficient fits with two 
parameter model (M35 and NGC 2548, see Paper II and 
Section 5.1.7). The maximum fraction of aligned stars 
favored by the three-parameter model is 0.3 (for NGC 
2548). For the cluster-by-cluster insights, if we treat 
M35 and NGC 2548 as possibly partially aligned, this 
yields a fraction of 2/10 or 20% consistent with partial 
alignment. 

'The small number of stars and potential for system- 
atic error contributing to these classifications highlights 
a key difference of our star-by-star discussion in Section 
5.4. With the star-by-star approach, these two clusters' 
stars are separated into isotropic and aligned subsets ac- 
cording to the f value of the three-parameter model. As 
a result, the influence of the small number of potentially 
aligned stars is reduced, since instead of a single clas- 
sification for each cluster, our results are weighted by 
the overall number of stars in the aligned and isotropic 
subsets. 


5.3.3. Isotropic or Moderately Aligned Clusters 


We have classified eight clusters having simi- 
lar peninsula-shaped posterior plots suggesting either 
isotropic or moderately aligned spins. For any single 
cluster, we cannot demonstrate a preference for one of 
these scenarios or the other. However, the ensemble 
of similar constraints facilitates further analysis. The 
possible explanations for this eight-cluster subsample 
are that their spin-axis distributions are all moderately 
aligned, all isotropic, or a mix of both. 

If all of the eight clusters in the peninsula-plot sub- 
sample have isotropic stellar spins, then at least 8096 of 
the 10 constrained clusters in our study have isotropic 
spins, while the remaining 2096 may be partially aligned 
(yet still with sizeable isotropic subsets, see Section 5.4). 
However, there remains the possibility that some or all 
of the isotropic/moderately aligned clusters are mod- 
erately aligned. We can put a limit on this number by 
comparing the ensemble of œ values associated with each 
cluster's moderate alignment to an expected distribution 
of a values. 
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To estimate the ensemble distribution of œ values re- 
sulting from the seven-cluster subsample, we used the 
PPDs from an MCMC run while limiting A < 65° to 
exclude the region of parameter space corresponding to 
isotropy or near-isotropy (for which « is not well con- 
strained, see Figure 5). We approximated these PPDs as 
normal distributions and mapped values « 0° or > 90? 
into the 0-90? range (e.g. —10? became 10?, and 100? 
became 80?). 'This "compressed" normal distribution 
well approximated the PPDs of a for each cluster, ac- 
curate to within a few degrees at the 5th, 16th, 32nd, 
50th, 84th, and 95th key percentiles of the distribution. 

We sampled o values from each approximated dis- 
tribution, choosing a number of values proportional to 
each cluster's maximum log-posterior value. We thus 
weighted sampling by our confidence in each cluster's 
result, since the clusters with the largest In(post) values 
received the greatest representation. Combining all the 
results of this sampling across clusters to create a sin- 
gle posterior distribution for a, we determined that the 
ensemble's o values in the aligned region of parameter 
space are > 32° at 95% confidence and range from ~ 37? 
to 73° at 80% confidence. 

Having estimated the range of observed «a values, our 
ensemble analysis continued by establishing an expected 
a distribution to compare to the data. We considered 
two scenarios: random orientations of all clusters’ mean 
inclinations, and uniform alignment. This latter sce- 
nario is simple but unlikely, because it requires a strong 
connection between isolated clusters amid the turbulent 
interstellar medium (e.g. Falceta-Goncalves et al. 2014). 
It would also be inconsistent with the spread of LOS 
angles to each cluster in our sample (Galactic latitudes 
between ~ —79° and 33° (standard deviation of 29°) and 
Galactic longitudes between ~ 15° and 274° (standard 
deviation 76°). 

If all clusters have their stellar spins oriented in the 
same three-dimensional direction, we would expect to 
observe a wide range of a values due to the varying 
LOS to each cluster. Instead, the comparable o values 
in our subsample contrast with the wide range of LOS 
directions, suggesting that not all the stellar spin-axis 
orientations in our subsample are intrinsically aligned 
with each other on a cluster-to-cluster scale. We there- 
fore make the assumption of random a values among 
any aligned clusters, with the cumulative probability of 
a given a value described by 


Pla) = 1 — cosa (3) 


(see Jackson & Jeffries 2010). 
We continued by assuming that all clusters’ œ values 
were within the 80% confidence range between 37° and 


73°. We then took the difference of cumulative proba- 
bilities from Equation 3 at both of these values to find 
that the probability of a given aligned cluster having an 
o value within this range is ~ 0.5 (assuming random 
orientations between clusters). The independence be- 
tween each cluster’s spin-axis orientation further yields 
the probability of n clusters having a values between 37° 
and 73° to be 0.5” from a simplified binomial probability 
distribution. 

For no fewer than 4 clusters to have o values within 
this range (with the remaining 4 isotropic, where a be- 
comes meaningless), the probability is p ~ 0.06. Thus, 
the probability of no more than 4 clusters being aligned 
in this a range is 1 — p ~ 0.94. Since the o range we 
used to determine this probability has a confidence level 
of 80%, we place a (0.8*0.94) = 75% confidence limit on 
the number of clusters having moderately aligned distri- 
butions at no more than four. It follows that no fewer 
than four of the eight clusters in this subsample have 
isotropic spins at 75% confidence. 


5.3.4. Constraints on Spin Alignment 


As a result of the above discussion, we make the fol- 
lowing constraints on spin alignment on a cluster-by- 
cluster scale at 75% confidence. Since at least four 
our of eight isotropic/moderately aligned clusters have 
isotropic spins, combining them with the two potential 
partially aligned clusters places a lower limit of four out 
of ten clusters (40% having isotropic spin-axis orien- 
tations. We also place an upper limit on spin-aligned 
clusters at six out of ten (60%). Of these six, up to 
four may have their entire spin-axis distribution moder- 
ately aligned, while up to two others are modeled as 
partially aligned and paired with an isotropic subset 
of stars. We can only place lower limits on the num- 
ber of clusters with isotropic spin-axis orientations (and 
upper limits on the number of totally/partially aligned 
clusters) because it remains possible that the latter are 
isotropic and either masked by the degenerate parame- 
ters of their peninsula-shaped posterior plots or experi- 
encing systematic error in their sin į values that is overfit 
by the three-parameter model. The number of isotropic 
clusters increases to 80% if the entire eight-cluster sub- 
sample has randomly oriented spins. 

For the eight-cluster subsample, we also performed an 
ensemble analysis of the spread half-angle posteriors in 
an analogous way to the mean inclinations. This time, 
we modeled the posteriors as normal distributions each 
having a mean A = 90°, the approximate PPD peak for 
all nine clusters. As for the ensemble mean inclination 
distribution, for values of each normal distribution < 0° 
or > 90°, we mapped them such that they were within 
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the range of 0-90°. This approximation matched the key 
percentile values of each A PPD to within a few degrees. 
From the combination of posterior-weighted draws from 
each of these distributions, we determined that > 25° 
at 95% confidence for eight-cluster subsample. 

While these constraints rule out tight alignment 
among the eight isotropic/moderately aligned clusters, 
the best-fitting three-parameter models for M35 and 
NGC 2548 favor tight alignment. However, the parame- 
ters are poorly constrained compared to the eight-cluster 
subsample due to the small number of stars (~ 11, see 
Section 5.4) modeled as aligned. In the next section, 
we discuss the interpretation of these results on a star- 
by-star basis, which can provide results less biased by 
potentially outlying sin ? values on each cluster’s classi- 
fication. 


5.4. Combined Insights: Star-by-Star 


Overall in the 15-cluster sample, we determined 404 
projected inclinations. The tally by cluster is in Ta- 
ble 5. The five unconstrained clusters contain 65 stars 
(~ 16% of the sample). This leaves 339 stars which 
are part of clusters with constrained inclination distri- 
butions. The eight-cluster isotropic/moderately aligned 
subsample contains 248 stars out of those 339, or 73%. 
At the lower limit of four out of eight isotropic spin-axis 
distributions in this subsample, we estimate that at least 
~ 124 stars are part of isotropic distributions. Since we 
cannot determine which clusters of the nine may have 
spin alignment, this value is not exact, but variable to 
a degree proportional to the number of sin? values for 
each cluster. The remaining ^ 124 stars may be part of 
an aligned distribution, or may also have isotropic spins. 

'The most significant difference using the star-by-star 
approach compared to cluster-by-cluster is the treat- 
ment of clusters modeled as partially aligned. While 
these clusters have a single classification determined by 
their disagreement with the two-parameter model, the 
star-by-star approach can consider the fraction of stars 
aligned versus those modeled as isotropic. For M35 and 
NGC 2548, best-fitting aligned fractions of f ~ 0.06 
and 0.3 mean that a majority of these clusters! stars are 
modeled as isotropic. When we break down these stars 
into aligned and isotropic categories, we find that ~ 11 
stars out of the entire 339-star constrained sample are 
modeled as aligned. The remaining 80 stars in these 
two clusters are part of the isotropic background in the 
three-parameter model. 

This result means that overall, 97% of stars with sin i 
determinations in our sample are consistent with being 
part of isotropic distributions. It is important to note 
that a particular degree of moderate alignment would 


present the same background sin? CDF as isotropy for 
the partially aligned clusters, a degeneracy analogous 
to the two-parameter model. The remaining 11 out of 
339 stars (~ 396) are modeled as aligned, although the 
potential for systematic error discussed above for these 
clusters and the small sample size of these stars prevents 
interpreting that they are indeed physically aligned. 


5.5. Implications for Star Formation 


The consistency of our results with a majority of 
isotropic spin-axis orientations suggests that among our 
sample of clusters with masses ~ 200-4000 Mo and 
ages ~ 5-700 Myr, turbulence dominates ordered rota- 
tion in the clumps from which protostellar cores form. 
For isotropic distributions simulated by Corsaro et al. 
(2017), there is a negligible amount of kinetic energy 
in rotation compared to turbulence. The remaining 
possibility of moderately aligned distributions leaves 
open the possibility that rotational kinetic energy for 
some clumps is > 10% of their turbulent kinetic energy. 
Nonetheless, turbulence dominates rotation for at least 
some of the clusters in our sample. This is in agreement 
with recent measurements of a star-forming filament's 
specific angular momentum that found a value compa- 
rable with protostellar cores in other regions, supporting 
a star's inheritance of its spin-axis orientation from the 
surrounding turbulence (Hsieh et al. 2021). 

'The mass of a cluster, if proportional to the mass of 
its progenitor clump, may have an important influence 
on the alignment of stellar spins. A clump with greater 
mass has a larger amount of angular momentum to im- 
part onto individual protostars, potentially countering 
the randomization of turbulence with greater efficiency 
than in less massive clusters. To present our results in 
light of this possibility, we plot in Figure 13 the spread 
half-angle (A) posteriors we determined for each of the 
ten constrained clusters versus the cluster mass and un- 
certainties from the literature. We plot M35 and NGC 
2548 using red squares to highlight their inconsistency 
with isotropic spins, as opposed to the rest of the sample 
labeled by black circles. 

'The minimum mass among the two clusters inconsis- 
tent with isotropy is ~ 800Mo, while the overall sam- 
ple contains clusters with masses as small as ~ 200Mo. 
A dependence of spin alignment on cluster mass could 
manifest as a trend in a diagram like Figure 13. How- 
ever, owing to the uncertainties and sample size of clus- 
ters in our study, we cannot identify any significant mass 
dependence of spin alignment. Future work may add 
value to this plot increasing the number of points and 
reducing their uncertainties. The reported strong spin 
alignment in the massive open clusters NGC 6791 and 
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NGC 6819 (Corsaro et al. 2017), with masses ~ 5000 
and 2600 Mo, respectively (Platais et al. 2011; Kalirai 
et al. 2001), further motivates the consideration of open 
cluster mass when studying the spins of their stars. 

While Corsaro et al. (2017) did not model aligned sub- 
sets of stars in a cluster, observations and simulations 
suggest that star formation involves hierarchical sub- 
structure even within a single clump (e.g. Judrez et al. 
2019; Lim et al. 2021; Torniamenti et al. 2022). The 
presence of such substructure during star formation sug- 
gests that the stars in a cluster may originate from many 
structures that were once distinguishable, and that spin 
alignment is more likely to manifest on these smaller 
scales compared to the overall cluster. In addressing 
clusters where a significantly populated subset of stars 
is discrepant with isotropic spins, the three-parameter 
model will continue to be useful in future studies. 

That the majority of our results are consistent with 
isotropy also supports the turbulence-driven misalign- 
ment of protostellar rotation axes with their magnetic 
fields. Numerical simulations by Kuznetsova et al. 
(2020) found that magnetic fields tend to be aligned 
among protostars. If their rotation axes are isotropic, 
however, then these protostars will tend to be magnet- 
ically misaligned and will experience a reduced level of 
magnetic braking. This misalignment reduces the effect 
of the magnetic braking “catastrophe,” in which a level 
of braking strong enough to decrease a protostar’s spe- 
cific angular momentum to observed levels also prevents 
the formation of protoplanetary disks. The isotropy- 
induced prevention of excess angular momentum diffu- 
sion supports the formation of protoplanetary disks on 
a scale similar to their observed properties (Gray et al. 
2018; Zhao et al. 2020). 


5.6. Implications for Radius Inflation 


The systematic offsets observed in some fits to the 
data in Figure 6 suggest the possible presence of radius 
inflation to explain our model’s underestimated sin į val- 
ues (see Section 3.5). By attributing this offset to ra- 
dius inflation and assuming isotropic spins, we can use 
the SPOTS models to estimate the amount of radius 
inflation in each cluster. 

We assumed isotropic spins in the 11 clusters for which 
isotropy is a possible solution (the eight with peninsula- 
shaped likelihood plots along with the unconstrained 
ASCC 16, Pozzo 1, and BH 56). We then shifted the 
isotropic model for each of these clusters towards larger 
sini values (rightward on the plots) by the mean frac- 
tional radius inflation amount predicted by the cluster’s 
SPOTS isochrone for all fractional spot coverage values. 
The minimum x2, yielded by adopting these values of 
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Figure 13. Spread half-angle posteriors versus cluster mass 
for the ten constrained clusters in our sample. Black cir- 
cles represent the eight clusters consistent with isotropic 
spins, while red squares correspond to NGC 2548 and M35, 
whose spin-axis distributions are modeled as inconsistent 
with isotropy. 


fspot revealed the amount of radius inflation necessary 
to mitigate the visible systematic offset between data 
and model for each of the 11 clusters. 

Figure 14 plots the mean fractional radius inflation 
amount corresponding to each cluster’s optimal isotropic 
fit as a function of its age in Myr. The points are color- 
coded by the fractional spot coverage value of the opti- 
mal result, and the error bars come from the standard 
deviation of the SPOTS versus MIST radius predictions. 
We find that, under our assumption of isotropic spins, 
the typical radius inflation among our sample’s mostly 
Sun-like stars on slow-rotator sequences is ~ 10% for 
clusters younger than ~ 100 Myr and drops to no more 
than a few percent for older clusters. The spot covering 
fractions associated with this amount of radius inflation 
range from 0 (one cluster best fit by the original MIST 
model) to 0.5 (two). 

Our radius inflation result is quantitatively similar to 
that of Somers & Pinsonneault (2015), who estimated 
up to 10% radius inflation in pre-main sequence stars 
and up to 9% for partially convective stars on the zero- 
age main sequence. The age threshold we determine for 
~ 10% radius inflation is also in agreement with this 
work, since their analysis of the ~ 100 Myr-old Pleiades 
found that the cluster’s stars on the slow-rotator se- 
quence were consistent with unspotted models yielding 
no radius inflation. Finally, Section 2.2.2 of Somers & 
Pinsonneault describes the observational evidence that 
pre-main sequence stars can have spot covering fractions 
up to ~ 50%, supporting the maximum fraction that we 
estimate in our analysis. 
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Figure 14. Mean (with standard deviation error bars) of the 
fractional radius difference between MIST isochrones and the 
SPOTS isochrone that achieves best fit to sini data assuming 
an isotropic distribution. All isotropic/moderately aligned 
and two of three unconstrained clusters (ASCC 16 and BH 
56) are plotted. Points are color-coded by the fractional 
starspot coverage associated with the best fit. 


5.7. Simulations Informing Future Studies 


To inform future studies searching for aligned subsets 
of stars in clusters, we performed Monte Carlo simu- 
lations of many partially aligned inclination distribu- 
tions. Our goal was to estimate the ability of the three- 
parameter model to confidently identify different frac- 
tions of aligned stars in clusters. We generated sin i val- 
ues using the three-parameter model, assuming random 
uncertainties corresponding to the mean values of the 
overall ensemble of measurements in this study (~ 11%, 
5% and 1% uncertainties respectively in vsini, P and 
R). In our simulations, the possible values of f ranged 
from 0.05-1.0 in increments of 0.05. For each value of f, 
we drew 5000 mean inclination values from a distribu- 
tion of isotropic Qaniso values, once again assuming no 
correlated spin alignment between clusters. 

We fixed the spread half-angle Aaniso in each simu- 
lation to a value informed by magnetohydrodynamic 
simulations of star formation. Considering that our 
work does not find evidence for the alignment of pro- 
tostellar rotation axes and magnetic fields, we consulted 
these simulations to estimate the minimum level of mis- 
alignment that permits the formation of protoplanetary 
disks. Based on the results of Hennebelle & Ciardi 
(2009), a misalignment of ~ 10? can moderate mag- 
netic braking at a protostellar mass-to-flux ratio con- 
sistent with observations to enable disk formation. We 
adopt twice this value, Aaniso = 20°, to represent the 
cone within which stellar spin axes are uniformly dis- 
tributed in our simulations. In doing so, we estimate a 


lower limit for the degree of spin alignment that a subset 
of stars can display. 

We fit each set of sini values with two CDFs: one 
for an isotropic distribution, and one for the three- 
parameter values used to generate each inclination dis- 
tribution. We compared the ratio of posterior values de- 
scribing each CDF's fit to the simulated data. A poste- 
rior value 19 times greater for the three-parameter model 
CDF than the isotropic CDF corresponds to a 9596 con- 
fidence determination of partial alignment over isotropy. 
We performed this process for simulated sample sizes of 
Nein i = 10, 30, and 100. 

Figure 15 plots the results of these simulations. The 
left panel plots the mean confidence level favoring par- 
tial alignment over isotropy as a function of the frac- 
tion of aligned stars f. The right panel plots the frac- 
tion of clusters Felusters for which partial alignment was 
favored over isotropy at > 95% confidence. Since ev- 
ery cluster in these simulations was generated from the 
three-parameter model with f > 0, the isotropic model 
is never the correct solution. The degree to which it is 
comparable with the appropriate three-parameter model 
solution offers insight into the outcomes of future studies 
of stellar spins in open clusters. 

We found that the minimum value of f required to 
produce a 95% confidence determination of partial align- 
ment, averaged over Qaniso values, was 0.9 for a distribu- 
tion of 10 sinz values and 0.85 for both 30 and 100 values 
(see colored symbols in the left panel of Figure 15). This 
suggests that a large fraction of stars needs to be aligned 
in a randomly selected cluster to yield a strong signa- 
ture of partial alignment. We also determined that the 
fraction of clusters with 95% confidence identifications 
of partial alignment is zero until f = 0.15, after which 
it rises to ~ 0.1 at f = 0.5 and increases until approach- 
ing ~ 0.25 as f approaches 1 (see colored symbols in the 
right panel of Figure 15). 

Approximately 9896 of the partially aligned spin-axis 
distributions detected at > 9596 confidence had Aaniso < 
45°. This shows that aligned spins with orientations 
closer to edge-on will be very difficult to confidently de- 
tect given the assumed alignment model, measurement 
uncertainties and rotation velocity distribution. For a 
collection of ten clusters each having half of their stars' 
spins partially aligned (with the other half isotropic), 
our simulations suggest that one will be identified as 
partially aligned by the three-parameter model at 9596 
confidence versus isotropy. Improving on the current 
~ 10% precision for vsini values, which was the dom- 
inant source of uncertainty in our study compared to 
P and R, would permit a greater number of high- 
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confidence identifications of partial alignment in a given 
cluster sample. 

We found little dependence of the mean confidence 
level and Fetusters on the sample size (10 € Nging < 100). 
There is a 5% decrease in the minimum f to achieve a 
mean confidence of 95% between 10 and 30 sini values. 
In practice, the possibility for apparent alignment due 
to a small sample of stars and systematic error effects 
on the spectro-photometric parameters informing sin? 
determinations supports the analysis of as many homo- 
geneously obtained data points as possible. However, 
returns diminish between ~ 30 and 100 sin i values, and 
upon reaching this level observing time becomes better 
spent studying the stars of another cluster. 


5.8. Kinematics 


The mean tangential velocities for ASCC 16, Alpha 
Per, and NGC 2547 are discrepant with zero at greater 
than 95% confidence. In each case, the direction of the 
motion is counterclockwise on the sky. Since the Galac- 
tic latitude of each of these clusters is negative, their 
counterclockwise rotation opposes the galactic disk an- 
gular momentum. This is consistent with models pre- 
sented by Kroupa et al. (2022). Since the inclination 
distributions of these clusters are either unconstrained 
or potentially isotropic, we cannot definitively compare 
the mean inclination of their spin-axis distributions to 
the orientation of their bulk rotation. 

Six clusters display significant radial kinematics. Of 
these, five (Collinder 69, ASCC 16, ASCC 19, Gulliver 
6, and Blanco 1) are directed outward from cluster cen- 
ter. Collinder 69, the youngest cluster in our sample 
(~ 5 Myr) also shows a significant signature of rotation 
about an axis perpendicular to the LOS, and it has a 
significant discrepancy between velocity dispersions in 
the radial and tangential directions. The expansion of 
young clusters is expected as a result of the expulsion 
of their gas (e.g. Kroupa et al. 2001), and Kuhn et al. 
(2019) found signatures of expansion in 75% of a sample 
of 28 clusters and associations between 1-5 Myr in age. 

One cluster (Alpha Per) has its radial motion directed 
inward. Makarov (2006) also found a signature of con- 
traction for the stars in this cluster by comparing their 
LOS and astrometric velocity dispersions. As we de- 
tailed in Section 5.8 of Paper II, the mean radial mo- 
tions we report should be interpreted with caution due 
to the potential bias introduced by the G < 18 magni- 
tude cutoff for cluster members in our sample. 


6. CONCLUSIONS 


We have analyzed stellar rotation periods, radii, and 
projected rotation velocities to determine the orienta- 
tion of stellar spin axes in 11 open clusters. Using two- 


and three-parameter models imagining stellar spin axes 
as uniformly distributed in a cone of varying mean in- 
clinations and spread half-angles, we place constraints 
on the inclination distributions of six of these clusters. 
We combine these results with those of four additional 
clusters from Papers I and II, resulting in a combined 
sample of 15 clusters. We perform an ensemble analysis 
on ten of these clusters with constrained spin-axis distri- 
butions (and their associated 339 projected inclination 
values) to provide insight into the initial conditions of 
star formation. 

We find that eight clusters (80%) are consistent with 
isotropic spins. These clusters have a mean mass ~ 1100 
Mo (range ~ 200-4000 Mo) and mean age ~ 100 Myr 
(range ~ 5-700 Myr). Among these clusters, this re- 
sult suggests the dominance of turbulence over ordered 
rotation in their progenitor clumps of molecular gas. 
Isotropic rotation axes also support a resolution to the 
protostellar magnetic braking “catastrophe” via reduced 
angular momentum diffusion by a misaligned magnetic 
field. 

An isotropic model degeneracy with parameters de- 
scribing moderate alignment prevents us from defini- 
tively reporting these clusters to have isotropic spins. 
Instead, we motivate and use the assumption of uncorre- 
lated spin-axis orientations between clusters to strongly 
rule out that all eight clusters are aligned and to estab- 
lish an upper limit of four potentially aligned clusters at 
75% confidence. This limit corresponds to no fewer than 
four isotropic clusters (50%) in the overall ten-cluster 
sample. Using the two-parameter cone model, we con- 
strain the mean inclination of potential spin alignment 
in Collinder 69, ASCC 19, NGC 2547, Alpha Per, Blanco 
1, NGC 2516, the Pleiades, and Praesepe to be between 
37° and 73° at 80% confidence and > 32° at 95% con- 
fidence. We also constrain the spread half-angle to be 
> 25° at 95% confidence. For these eight clusters, we 
thus rule out tight alignment and mean inclinations ori- 
ented along the LOS. 

We find that two clusters’ inclination distributions are 
better fit by a three-parameter model that superposes 
aligned and isotropic subsets of stars. However, the com- 
bination of the small number of stars modeled as aligned 
and the potential for systematic error in their sin ¿į values 
prevents an interpretation of spin alignment. The ma- 
jority of these stars’ fellow cluster members are modeled 
as isotropic, resulting in a 97% fraction of stars consis- 
tent with isotropic spins across our ten-cluster sample. 

Future studies of stellar inclinations in open clusters 
could increase the number of clusters studied and fur- 
ther reduce systematic error in order to improve con- 
straints on spin-axis distributions. Having a greater 
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Figure 15. Partial alignment simulation results for sample sizes Nsini = 10, 30, and 100. Left panel: mean confidence in 
partial alignment over isotropy versus fraction of aligned stars f. Colored symbols indicate mean confidence > 95%. Right 
panel: fraction of clusters Folusters identified as partially aligned (at 95% confidence over isotropic) versus f. Colored symbols 


indicate a fraction > 0. 


number of clusters can strengthen the lower limit on 
isotropic versus moderately aligned clusters. A greater 
number of stars per cluster can be beneficial as well, al- 
though with diminishing returns on spin-axis constraints 
beyond sample sizes larger than ~ 30 stars (e.g. Jack- 
son & Jeffries 2010). Therefore, increasing the number 
of clusters is more useful to improve the lower limit of 
isotropic clusters. A larger cluster sample size will also 
facilitate inquiry into the mass dependence of spin align- 
ment. 

We recommend multi-object spectrograph studies tar- 
geting Gaia-based open cluster members across the sky. 
The projected rotation velocity remains the least avail- 
able and most inhomogeneous quantity of the three pa- 
rameters of Eq. 1. The above proposed studies would 
relieve that shortage and provide more homogeneous 
vsini results. In addition, further understanding mag- 
netic radius inflation (see Section 5.6) would reduce a 
key systematic error that we encountered, especially for 
the youngest clusters in our sample. Increasing the num- 
ber of clusters shares priority with reducing systematic 
error in all parameters of Eq. 1 for future work. 

As illustrated by the 11 stars modeled as aligned in 
M35 and NGC 2548, even a small number of stars can 
lead to a statistically significant departure from isotropy, 
which increases with the fraction of aligned stars (see 
Section 5.7). Increasing the number of sini values for 
clusters that are not well fit by the two-parameter model 
will reduce the possibility of illusory alignment and yield 
tighter constraints on physical alignment (again facing 
diminishing returns beyond ~ 30 stars). The orienta- 


tion of partially aligned subsets has a large influence 
on their detectability, with pole-on orientations much 
easier to differentiate from an isotropic subset. A key 
improvement for future studies will come from increas- 
ing data quality and minimizing systematic errors, help- 
ing to confidently determine whether signatures of align- 
ment are physical or instrumental in nature. 

While TESS rotation periods were plentiful and es- 
sential to this work, the upcoming PLATO mission's 
smaller angular resolution will reduce the ambiguity due 
to blending with which we contended (Rauer et al. 2014). 
PLATO will also enable asteroseismic studies of stellar 
inclination that, while subject to their own unique bi- 
ases, provide additional direct measurements to compare 
with the spectro-photometric inferences of sin i. 
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